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In a system with an even number of SU(2) spins, there is an overcomplete set of states — consisting 
of all possible pairings of the spins into valence bonds — that spans the 5 = Hilbert subspace. 
Operator expectation values in this basis are related to the properties of the closed loops that 
are formed by the overlap of valence bond states. We construct a generating function for spin 
correlation functions of arbitrary order and show that all nonvanishing contributions arise from 
configurations that are topologically irreducible. We derive explicit formulas for the correlation 
functions at second, fourth, and sixth order. We then extend the valence bond basis to include 
triplet bonds and discuss how to compute properties that are related to operators acting outside the 
singlet sector. These results are relevant to analytical calculations and to numerical valence bond 
simulations using quantum Monte Carlo, variational wavefunctions, or exact diagonalization. 



I. INTRODUCTION 

The traditional way to represent the quantum states 
of a system of S = 1/2 spins is to introduce a basis of 
S z eigenstates. Each state corresponds to a particular 
assignment of "up" or "down" to each spin in the lattice. 
This basis is complete and orthonormal. Another useful 
basis, dating back to Rumer and Pauling in the 1930sji*£ 
is the so-called valence bond basis, in which the states of 
the system are represented by partitions of the spins into 
pairs forming singlets. A valence bond (VB) state, 
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FIG. 1: The low-energy singlet sector of quantum antiferro- 
magnets can be described in either the S* or valence bond 
basis. 



is a product of singlets, 
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in which each site label appears only once. Pictorially, 
a VB state can be represented as a system of hard-core 
dimers, where each lattice site (spin) belongs to exactly 
one dimer, as illustrated in Fig. ^ The VB basis spans 
the singlet sector of the Hilbert space, but is massively 
overcomplete. It is also highly nonorthogonal, having the 
unusual property that every two states in the basis have 
nonzero overlap. 

The first important work in this basis was the calcu- 
lation of the ground-state energy per spin of the infinite 
quantum Heisenberg chain by Hulthen 3 (building on the 
earlier work of Bethel). Hulthen also determined the 
eigenvalues and eigenstates of small finite chains of up to 
ten sites using the subset of "noncrossing"«i valence bond 
states (which form a complete basis). Majumdar and 
Ghosh later performed a similar calculation for the one- 
dimensional spin chain with nearest- and next-nearest- 
neighbour interactions^ 

The valence bond basis has a special connection to 
the paramagnetic states of interacting antiferromagnets. 
Fazekas and Anderson introduced the resonating valence 
bond (RVB) picture to describe a possible spin liquid 
in frustrated magnetspi Anderson later proposed that a 



doped RVB state may describe the cuprate supercon- 
ductorsi^ This suggestion spurred wide interest in the 
VB description of antiferromagnets and possible exotic 
ground states that are naturally described in terms of 
VB states. In particular, VB states have proved useful 
as variational states, where the variational freedom lies 
in the distribution of valence bond Iengths2ii2iii or in a 
set of projected-BCS coefficientsi 12 i 13 i 14 

Since valence bond states have total spin invariance 
built in, they are also a natural choice for exact diag- 
onalization within the low-energy S = matrix block. 
Such calculations have been used to study RVB states 
with a restriction of the VB basis to include only a sub- 
set of states with short dimers which should be a 
good approximation for systems with only short-range 
spin correlations. A generalization of this procedure — in 
which the valence bonds are antisymmetrized products of 
individual atomic eigenstates — has been widely adopted 
by computational chemists for use in molecular quantum 
mechanics ]F4§ 

Variational calculations suggest that the collinear Neel 
ground state of the d-dimensional (d > 1) Heisenberg 
antifcrromagnet is well-described by a superposition of 
VB configurations whose distribution of bond lengths 
exhibits l/r p powerlaw behaviour, with p < 5 in two 
dimensions £ In two dimensions, it is also believed that 
there are magnetically disordered states with energy very 
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close to that of the ordered ground state: 9 Thus, compet- 
ing interactions may favour RVB states with short-ranged 
bonds and no long-range magnetic order i 1 ?' 20 They may 
also favour states with valence bond solid (VBS) or- 
derfSliSSi in which the translational symmetry is broken 
but the spin-rotational symmetry remains intact. To 
date, the only firm confirmation of an RVB state is in 
simplified quantum dimer models, 23 in which only the 
dimer degrees of freedom are retained and the spin de- 
grees of freedom associated with the dimers are neglected. 

In most cases, possible exotic spin states^2^SS*21 sim- 
ply cannot be studied in detail with unbiased methods 
such as quantum Monte Carlo (QMC), because of the 
infamous negative-sign problem. QMC simulations of 
quantum spin systems have traditionally been carried out 
in the standard S z basis. Although the feasibility of a 
Monte Carlo projection for improving a variational VB 
state was demonstrated more than fifteen years agOr& 
the use of the VB basis in QMC studies has been very 
limited so fariS As it turns out, neither the overcomplete- 
ness nor the nonorthogonality of the valence bond basis is 
an impediment to carrying out QMC simulations in prin- 
ciple. In fact, even without a good variational state as a 
trial state, the ground state can be completely projected 
out starting with, e.g., an arbitrarily chosen basis stated 
using importance sampling and a simple local updating 
scheme. This method delivers performance comparable 
to the current state of the art for 5 z -basis QMC calcula- 
tions (i.e., stochastic series expansion 31 and world line** 
methods with loop-cluster updates) . It was also noted in 
Ref. |3(] that the VB projector method expands the class 
of models that are sign-problem-free. A host of isotropic 
SU(2) invariant models with multi-spin interactions (in- 
volving four, six, etc. spins) can be studied, which in 
spite of not being frustrated in the standard sense (an 
odd number of antiferromagnetic interactions around a 
closed loop on the lattice) do give rise to a sign problem 
in standard methods. Thus, VB simulations open new 
opportunities to study ground state phases and transi- 
tions in quantum spin systems. 

In addition to these new opportunities for QMC studies 
in the VB basis, we also believe that there is more to do 
variationally. There was not much follow-up on the pio- 
neering calculations by Liang, Doucot, and Anderson, 9 
and with the increase in computer performance since 
that time, it is now feasible to consider more complex 
wave-function optimizations. 11 Although frustrated sys- 
tems also cause sign problems in variational calculations, 
it may still be possible to extract useful information from 
them. It may also prove fruitful to study variational VB 
wave functions for nonfrustrated multi-spin interactions. 

In both QMC and variational calculations, one would 
like to study a wide range of physical observables to char- 
acterize the ground state. Overlaps and matrix elements 
between two VB states \v) and \v') can be related to 
the structure of the closed loops that are formed when 
their corresponding dimer configurations are superim- 
posed. There are two standard results: (i) For a system 



of 2N spins, (v\v') = ±2 N o~ N where N < N is the 
number of loops. (v\v') is unity when the states \v) and 
\v') are identical and Nq is a maximum; its magnitude 
is halved each time a bond mismatch reduces the loop 
count by oneM^l (ii) (v\Si ■ Sj\v') = ±j{v\v') if i and 
j belong to the same loop and vanishes otherwise^ (In 
each case, the overall sign depends on the convention for 
assigning directions to the bonds.) We are not aware of 
expressions for more complicated matrix elements, e.g., 
those involving more than two spin operators or indi- 
vidual components of the spins. Nor are we aware of 
expressions for quantities that are associated with triplet 
excitations or for quantities, such as the spin stiffness, 
that cannot usually be written in terms of an equal-time 
correlation function. 

The goal of this paper is to provide a formal framework 
for organizing calculations in the valence bond basis and 
to present several new formulas relating the loop struc- 
ture of overlapping valence bond states to a wider range 
of physical properties of the system. We do this in a for- 
mal way, using bond operators equivalent to those first 
derived by Sachdev and Bhatt^ which in their usual 
context are associated with a fixed dimer configuration. 
Here, we generalize the creation and annihilation opera- 
tors that describe the states of the two-spin system to the 
many-spin case, allowing the operators to act between 
arbitrary pairs of spins. With one additional anticom- 
mutator rule, we find that we can use these operators to 
organize calculations in the overcomplete valence bond 
basis. As a result, we are then able to address higher- 
order spin correlations, such as (Si • Sj)(Sfc • S;), and 
demonstrate that their matrix elements are topological 
in nature, depending not only on the loop membership 
of the various site labels (as is the case for Si • Sj) but 
also on the overall connectivity of the loops with respect 
to the and (k,l) vertices. We introduce a gener- 

ating function whose derivatives produce a related class 
of cumulant function. A diagrammatic expansion eluci- 
dates the structure of these functions and greatly simpli- 
fies their calculation. We extend the valence bond basis 
to include triplet as well as singlet bonds in an effort 
to access the full Hilbert space. This allows us to com- 
pute various spin correlation functions component- wise, 
including two-spin operators of the form SfSJ and four- 
spin operators such as SfS^S^Sf and SfS^SlSf. We 
also derive an expression for the singlet-triplet gap and 
the spin stiffness. These results should be useful in the 
context of QMC, variational calculations, and in exact 
diagonalization. Our formalism may also find use in ap- 
proximate analytical calculations. 

The organization of the paper is as follows. In Sec. [H] 
we introduce the valence bond operator formalism that is 
used throughout. The three subsequent sections review 
how to construct 5 = valence bond states (Sec. IIIII) . 
how to describe their evolution under action by a hamil- 
tonian (Sec. liVfl . and how to compute their overlaps 
(Sec. 0). In Sec. IVII we begin to derive formulas for 
the isotropic spin correlation functions in a somewhat 



3 



naive way. We reproduce these results in Sec. IVIl| using 
a more sophisticated diagrammatic approach and then 
carry the calculation to higher order. In Sec. IVIIII we 
develop rules for valence bond states with triplet excita- 
tions; these are used in Sec. IIXI to compute correlation 
functions of components of the staggered magnetization. 
In Sec. |3 we discuss the Neel state and how it can be em- 
ployed as a reference state to compute the singlet-triplet 
gap fSec. IXIjl and the spin stiffness (Sec. IXIl|) . 



II. VALENCE BOND OPERATOR FORMALISM 



Two SU(2) spins labelled i and j have a total spin 
S = Sj + Sj satisfying ^S 2 = | +Sj • Sj. Thus, the usual 
Heisenberg interaction can be written as 



Hi 



= Si -Si- - = -S 2 - 1. 



(3) 



The eigenstates of Eq. (0 are states of well-defined total 
spin, having eigenvalues ^S(S + 1) — 1. We enumerate 
them as follows: one (S — 0) singlet |0) and three (S — 1) 
triplets |1>, |2), 1 3) with energy E» = -5»°. 

For bookkeeping purposes, it is helpful to introduce a 
boson representation for the spins, 



2 Yl h ™ CTss ' hl 



with 



(4) 



labels correspond to real physical directions in the stan- 
dard basis of R 3 . This is true in the sense that 5?|0) ~ \a) 
for a — I. 2,3. The key is that Eq. Q obeys err = it, 
which implies that SiX°||vac) = («'/2)%t j .|vac). 

We emphasize that the x% operators completely de- 
scribe the SU(2) <g) SU(2) ~ SO(4) degrees of freedom of 
the two-spin system. These operators obey the complete- 
ness relation 
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which follows from Ylu T s S ' T r'r = ^S sr 5 s ' r ' and the re- 
striction to one boson per site. The anticommutation 
relation 



[X? 7 ,xr/]|vac)=<nvac> 



(10) 



is inherited from the properties [&is,frj s /] = SijS ss > and 
6 ls |vac) = 0. 

Any valid (i.e., ^2 s b\ 8 bi s =1 preserving) operation on 
one or both of the spins has an SO (4) representation. By 
construction, the operator equivalence 
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Xij ' Xij 



holds, and elimination of x y i a Eq. @ yields 
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(note the single-occupancy constraint), and to define a 
valence bond operator 
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Of 
Xij Xij ■ 



(12) 



v 



(5) 



that creates eigenstates of ify out of the bosonic vac- 
cuum (i.e., = X^l vac ) an d Xijl vac ) = 0). The eigen- 
value equation Hij\fx) = E^\^l), now written as 

(6) 



^Xf> a c)=-^ x°j|vac}, 



determines the unknown coefficients r M . One possible 
solution to Eq. © is 



(t°,t) = (i<r 2 ,ia 3 ,±,— ia 1 ), 



(7) 



where = (1, er) denotes the four-vector consisting of 
the unit matrix and the three Pauli matrices; the result- 
ing states are 



I/') 



xf/l vac ) 



f^(ITUj) - Uiti» if/i = o 

r^(|TiTi> - Il4i)) if A» = 1 
-^(\UW + \idj)) if A* = 2 



(8) 



Other linear combinations of the triplet states are 
equally valid, but this choice has the advantage that the 



For an arbitrary bilinear operator O — Ylnv Xij Xij > 
two applications of Eq. l|l(Jf) will coax out the coefficient 
matrix, 0"" = (vac|[x M , [6, x vi ]] |vac). 

In the case of the spin operators themselves, one finds 
(SfY u = \tv(r^a a T v ) and {S?y = § tr(T u a a *T^). 
Evaluation of the traces leads to 



^i ^j ^(.Xij Xij Xij X 



^i &j * e Xij Xij- 



(13) 
(14) 



Equations (|13fl and l|14|) turn out to be very useful, but 
we should not regard them as the fundamental operator 
equivalence rules. They seem to suggest that two-spin 
operations are always quartic in x% i which is clearly not 
true [cf. Eq. (JTTJ]. 

There is an alternative to computing the coefficient 
matrix directly. For operators that transform in a known 
way, we can induce the same transformation in x% an( l 
equate the corresponding terms. Consider a spin rotation 
S a Rab( dn j S b of e rat ji ans a b ou t the axis n (|n| = 
1). In the boson language, this rotation is equivalent 
to the unitary transformation b — > Ub with U(9n) = 
e (ie/2)n-a = lcos (^/2) + i n . crsm(9/2). Now suppose 

that we rotate Sj and Sj by two different angles about 
the same axis: putting b. L — > U(9ir\)bi and bj — > U(6jXi)bj 
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into Eq. J5J) , we find that the rotation is equivalent to the 
valence bond operator transformation 



Xij 

where 9%j = 



X°cos(%/2) + n-XySin(^/2), 
— 9j is the relative rotation angle. 



(15) 



R 



If n is directed along the 3 axis, the rotation matrix 
ab - \ tr U ] a a Ua b has the form 



R(9e 3 ) 



cos 9 sin 9 ' 
— siv\.9 cos 9 
0/ 



(16) 



TABLE I: Operator equivalence rules 



spin basis 


valence bond basis 


Si " Sj 

■ Sj — -j 

Si Sj 
Si -)- Sj 

S<"^j + &i "Sj" 
H^i ^j ~~ <^j ) 


l^Xij Xij ~4-'%ij ' ~X,ij 

v ot 

Xij Xij 
^{XijXij ~ XijXij) 

ixtj x Xv 

~Xij Xij ~^ Xij Xij 
Xij Xij ~l~ Xij Xij 



Described in terms of raising and lowering operators, the 
transformation amounts to 



III. 



5 = VALENCE BOND BASIS 



and in particular, 



s i s j 
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'S+S- +e w "S+Sj 
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Consequently, the isotropic Heisenberg interaction, ex- 
panded in powers of 9ij , behaves as 



For a system of many SU(2) spins, the structure of 
the Hilbert space can be determined from the product 
rule % ® S = (S - ® (S + where S denotes the 
(2S + l)-degenerate spin-S 1 state. Thus, 
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Sj • Sj — > S, • Sj -9i 3 ( S',' .V, — 



N 



4 °ij( S t S j 
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According to Eq. I)15[l. the same transformation ap- 
plied to the valence bond operators gives 
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?tv° +Y 3t Y 3 ). 
l 3 A-^3 azj ^- % 3 ' 



(20) 



Comparing Eqs. I|19|l and 120|) allows us to make the iden- 
tification 



S7S: 



Of 3 i 3f 

•y 1 y - - -r- V V ■ - 
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(21) 



We have so far confined our discussion to a system of 
two spins. Nonetheless, everything derived up to this 
point applies equally well to any two spins of a many- 
spin system. In that sense, the results summarized in 
Tableware valid generally: we simply assert that there 
are operators x% associated with every pair of sites in the 
lattice. The only remaining step is to determine what the 
appropriate anticommutator algebra is for operators with 
one site index in common. A straightforward calculation 
shows that 

K , X$Xu] l vac > = \ T ^ X# l vac > . ( 22 ) 



where T x ^ v 
repeated index A is implied 



i tr T A V AI T p ^T ,y and summation over the 



27V times 



S=0 



The number, Cg , of 5-blocks that appears in the 
Hilbert space of 2N spins is given in Table [H] The num- 
ber of states in the singlet sector, Cq — j^fiC 2 ^) = 
n\(n+i) \ ' re P resen ts a small fraction of the total number 
of states, Es=o( 2S + 1 ) C s N) = ^ ■ 



TABLE II: Values of the coefficient C s , , representing the 
number of spin-S blocks in the Hilbert space of N pairs of 
SU(2) spins. 



1/2 1 3/2 2 5/2 3 7/2 



1 

2 
5 

14 
42 
132 



1 

2 
5 

14 
42 
132 



1 

3 

9 
28 
90 
297 



1 

4 
14 
48 
165 



1 

5 

20 
75 
285 



1 

6 

27 
120 



1 

7 

35 
144 



9 

24 10 
34 



To construct a valence bond state in the singlet sector 
[as per Eq. we simply group the spins into N pairs 



and, starting from the bosonic vacuum, act with xl] 
each pair The number of such states is 



for 



T N = 



1 {2N\ (2N - 2 



N\ V 2 J V 



(27V)! 
2 N (Nl) ' 



(24) 



This dwarfs the number of independent states in the sin- 
glet sector for large N (T N > C { N) ). 

As a consequence of the overcompleteness, linear 
combinations of valence bond states are generally not 
unique — an ambiguity related to the fact that linear in- 
dependence fails even at the level of two singlets: 



ot ot , Of Of 

Xil Xjk ' Xij Xkl 



Of Of n 



(25) 



In the case of two spins (N — 1), there is a single 
tiling (Ti = 1) corresponding to the unique singlet state 

(Cq 1 ^ = !)• I 11 the case of four spins (N = 2), the number 
of tilings (T 2 = 3) exceeds the number of physical singlet 

states (C^ = 2) by one. It is useful to eliminate this 
superfluous state. For concreteness, let us suppose that 
the lattice of 2N spins is described as the union A U B of 
sites A = {ii,i 2 , ■ ■ • Jn} and B = {ji,j 2 , ■ ■ -Jn}- This 
is simply a labelling trick and does not depend on the 
lattice being bipartite. We restrict the basis to include 
only those valence bonds that connect A sites to B sites. 
This is possible since any unwanted AA and BB bonds 
can be replaced by two AB bonds via Eq. (|25|l . 

The AB valence bond basis is intimately related to Sn, 
the symmetric group of degree N. There is a one-to-one 
correspondence between the bond tilings and the permu- 
tations that are the elements of Sn- This correspondence 
holds because every bond tiling can be thought of as a 
rearrangement of the B labels. That is, for each P £ Sn, 
there is a state 



P) = P f |vac) = 



n 

n=l 



Ot 



vac) 



(26) 



It follows that there are Tn — N\ such states. For exam- 
ple, Fig. [3 shows the 3! = 6 AB bond configurations for 
N = 3. These states are indexed by the permutations 



P=(l)(2)(3) 
P= (1 2 3) 



P = (l)(2 3) 
P = (3 2 1) 



P=(l 2)(3) 
P=(l 3)(2), 



written here in cycle notation. 

The restriction to the AB basis is useful in that it gives 
us a systematic way to fix the phase of each state. This 
is an issue because the singlet is a directed bond: x% — 
—X% an d Xji = Xij- We adopt the convention that the 
canonical form of the operator x% has i £ A and j e B, 
as in Eq. I|26|) . As we shall see in Sect. all overlaps 
of AB valence bond states are positive definite, which is 
not true in the unrestricted basis. 

One disadvantage to the AB basis is that the identity 
operator in the S — subspace has a complicated form. 



N = 1 



N = 3 



N = 2 







FIG. 2: The set of all possible bond tilings, T N = (2N)\/2 N N\ 
in number, is shown for N = 1,2,3. The restricted set of 
AB tilings (blue) numbers fjv — N\, which still exceeds the 
number of true singlet states. A and B sites are shown as 
filled and open circles, respectively. 



In the unrestricted basis, it is diagonal, 

2 N _ 



lo = 



v)(v\ 



(27) 



with a normalization Tn/Cq N \ but the process of elim- 
inating all AA and BB bonds from Eq. 1)2 7|l introduces 
offdiagonal terms. For N — 2 and N = 3, the identity 
operators are 



and 



(Q\ 



( 4 -1 -1 -1\ 
-1 4 0-1-1 



(28) 



-10 4 
-1 -1 
-1 -1 
V-l 



-1 -1 
4 0-1 
4-1 

-1 -1 4/ 



(Q\. (29) 
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In these examples, the rows and columns of the coefficient 
matrices are arranged to match the ordering of the per- 
mutations in Fig.EI The entries are (2N)\/2 N ~ 1 (N\) 2 -I 
along the main diagonal and or — 1 elsewhere, which is 
a consequence of there being at most one AA/BB set to 
unravel. For N > 3, the situation is more complicated, 
and we do not know of a simple general expression for 
lo- In most practical applications, however, such an ex- 
pression is not needed. 

It is possible to construct basis sets that are still more 
restrictive (i.e., having fewer states than the AB basis). 
There is a hierarchy of linear dependence relations [sim- 
ilar to Eq. l(2*5)l] for groups of three, four, and higher 
numbers of singlet bonds, which can be used to further 
eliminate unwanted states. When carried out to its fullest 
extent, this elimination procedure leaves a set of valence 
bond configurations equal in number to the actual num- 
ber of S — states. Explicit construction of the ba- 
sis in this limit can be accomplished by arranging the 
lattice sites «i,ji,«2, J2, • • • )ijv»jjv on a ring and keep- 
ing only the AB tilings that produce no bond crossings. 
Note that even this minimal set contains long bonds con- 
necting sites that are macroscopically separated, which 
suggests that bonds on all length scales are required for 
completeness. 

IV. EVOLUTION OF VALENCE BOND STATES 

Now suppose that we have a hamiltonian of the form 
H = JijHij — ^ KijkiHijHki + • ■ • (30) 

ij ijkl 

with second-, forth-, and potentially higher-order spin 
interationsi^i Here, ify is defined as in Eq. (|12|l . In order 
to understand how an arbitrary state evolves under H, 
we need to know how a given valence bond state evolves 
under repeated applications of —Hij. 

Depending on the circumstances, — H^ is either a di- 
agonal operation that leaves the overall configuration 
unchanged or an off-diagonal operation that maps two 
bonds to their complementary tiling: 

(i-s < -s i )M = M, (3i) 

(±-S i -S j )[l,i}\j,k] = [i,j][k,l}. (32) 

This result, which is true in the unrestricted basis, has 
long been knownpS. It will be instructive to see how it 
emerges within our formalism and how it is modified to 
accommodate the restriction to AB bonds. 

Let us start with the state \P) = P^vac), where P' is 
the operator string defined in Eq. (|26J) . Then, 

-Hij\P) = xM^lvac) = X-j[x° ,^]|vac>. (33) 

There are two possibilities to consider in evaluating the 
anticommutator: (i) If there is already a bond connecting 



the active sites (i.e., P^ — ■ ■ ■ xlj ■ ■ ■), then only two 
operators play a role in the anticommutator. In this case, 
Eq. ((THf is the appropriate rule to apply, and Eq. lt33l 
simplifies to 

X%X%\P) = \P)- (34) 

(ii) If the active sites are each connected elsewhere (i.e., 

P f = ■ • ■ Xa ■ ■ ■ xl] ■■ ■ f° r some k G A, I G B), then 
three operators are involved. This necessitates the use of 
Eq. 1221) and leads to 

X°jXij \P) = \xll JPl ■ ■ ■ X°j ■ ■ ■ X°k\ ■ ■ ■ Xilj PN l vac )- (35) 

The new state on the right-hand side is itself a valence 
bond state (different from \P)), since it consists of N 
operators and none of its 2N indices are repeated. If 
i E A and j G B, then this state contains only AB bonds. 

A reconfiguration of the bonds manifests itself as 
change in the permutation indexing the state. If i = i n 
and j = j m are in opposite sublattices, then Eqs. I|34|) 
and H35f) are summarized by the compact update rule 

-KjJP)=(~) Pn ' m \(Pnm)P), (36) 

where the term (Pn m) inside the ket is a 2-cycle acting 
on P (the effect of which is to swap the labels jp n and 
j m ). Equation (|36|l makes clear that the hamiltonian is 
an identity operation whenever it acts on a preexisting 
bond; otherwise it induces a single transposition and a 
multiplicative factor i. Note that the transposition itself 
depends on P. 

If, however, i and j are in the same sublattice, then 
the XijXki operators in Eq. I)35|l have to be eliminated in 
favour of AB bonds. This leads to 

-KiJP) = \\ p )- Q) n ' m \{Pn Pm)P) (37) 
and 

-H jnjm \P) = \\P) ^Kn m)P). (38) 

These "frustrating" interactions create a linear super- 
position of states, rather than just a rearrangement of 
bonds. They also introduce bond configurations with 
negative weight, violating the conditions of the Marshall 
sign theorem^ Figure illustrates the frustrating and 
nonfrustrating bond flips that can occur. 

If there are no frustrating interactions in the hamilto- 
nian then the ground state can be written as a superpo- 
sition of valence bond states — J2p c p\P) wfth all 
cp > 0. In that special case, Monte Carlo simulation is 
sign-problem- free. 



7 



Jm 




3r, 




Figure 0] illustrates the construction. Each of these loops 
is even in length — constituting a chain with equal num- 
bers of links from P and Q. The cycle decomposition is 
related to the loop membership of the site indices by 



QP 



Pm <£4> i n and j m in same loop 



n ~ m <^> i n and i m in same loop (41) 
Pn ~ Pm <^> ] n and j m in same loop 







FIG. 3: AB valence bonds obey a simple set of update rules 
when acted on with an isotropic spin interaction term. Here 
the update rules for 1/4 — Si ■ Sj axe summarized. The red bar 
denotes the interaction between the spins at sites i and j. The 
site labels follow the notation in Eqs. 1361 and 137H . Of the 
terms that involve a reconfiguration of the bonds (dashed out- 
lines indicate the previous bond locations), one corresponds 
to an exchange of jp n and j m and the other jp n and jp m . 



V. OVERLAPS AND MATRIX ELEMENTS OF 
VALENCE BOND STATES 

It follows from Eq. (|26[) that the overlap of any two 
valence bond states is equal to the vacuum-state expec- 
tation value of a length-2iV operator string: 



(Q\P) = <vac|xL, QN 



Al 2 ,702 A«i 



X Y 0t 
AllJP 



Of 

1 Ai 2 , JP2 



JQl 

ot 



X 



(39) 



The creation and annihilation operators can be shuffled 
past one another so long as they share no indices in com- 
mon. The maximal such rearrangement leaves the string 
grouped into several linked chains of operators; these are 
related to the disjoint cycles of the composite permuta- 
tion QP (we use the notation Q = Q _1 to denote the 
inverse permutation of Q) and define a set of directed 
loops (with QP defining the proper order), e.g., 



3 pi 



l QPl ~^ JPQPl 
~ * i-QPQPl ~ 



Hqp)"! = h- (40) 



QP 

We use the notation x ~ y to indicate that x and y are 
in the same cycle of QP. All the loops taken together 
form a set of closed paths covering the lattice. Hence, 



TV 



(42) 



fc=i 



where rik is the number of loops of length 2k (or the 
number of fc-cycles in QP). 

The value of the overlap can be computed by way of a 
simple decimation scheme^4^4 Each loop can be elimi- 
nated by repeated application of Eqs. I|1(J|) and l|22|) . spe- 
cialized to the case of p = v = fx = 0: 



and 



f Of Ofi 



Of I 



)Xki 



(43) 



(44) 



For a loop of length 2k, k — 1 applications of Eq. H44|) 
removes 2k — 2 links and yields k— 1 factors of 1/2. The 
final two operators are removed with a single application 
of Eq. H43fl . leaving (vacjvac) = 1. Hence, 



N 



(Q\P)=U 



1=1 



1 



fe-1 



2 k 



2£f=i«* 
2£i=i*">* 



N 



(45) 



where we have used Eq. I|42|) and defined the total number 

of loops No = J2k=i nfc - 

Similar arguments can be used to compute the ma- 
trix elements of a spin-rotation- invariant operator O. In 
general, such an operator acting on a valence bond state 
produces a linear superposition of states with modified 
bond configurations: 



6\p) = J20p,p' 
p' 



IP') 



(46) 



(The situation is more complicated if O generates states 
outside the 5 = sector; this is discussed further in 
Sect. IV111I ) Since the basis is overcomplete, the coef- 
ficients Op.pi arc not unique and Op t p< ^ (P\<J\P'), 
as would be the case in an orthonormal basis. Equa- 
tions i|45|) and l|46|l imply that matrix elements are re- 
lated to changes in the loop number: 



(Q\o\p) 
(Q\p) 



= 2J° 

p> 



p,P'2 N 'o- 



Nq 



(47) 
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FIG. 4: The top-left and top-right panels illustrate valence 
bond states \P) and \Q) corresponding to the permutations 
P = (1 2)(3)(4)(5)--- and Q = (1)(2 4 5 3) • • • ; the sites 
ii,£a, . . . , is (filled circles) and ■■■■,35 (open circles) are 

numbered accordingly (in no particular order). The over- 
lap (Q\P) between these states is literally that: a superpo- 
sition of the two bond configurations. Since the end point 
of one bond is always the starting point of another, a set of 
closed bond paths is formed. These are in one-to-one corre- 
spondence with the disjoint cycles of QP. In this example, 
QP = (1 2 4 5 3) • • • has a cycle of length 5 and a correspond- 
ing bond loop of length 10. 



Here, Nq is the number of cycles in QP and NU the 
number in QP'. 

Evaluation of Eq. I|47|l requires some general rules for 
how the number of cycles changes as P — ► P' . Since the 
permutation P' P can always be decomposed into a prod- 
uct of transpositions, it suffices to consider the situation 
where the two states differ by a single transposition. In 
that case, NL — iV^ = ±1 (see the Appendix^), since 
a transposition (n m) either merges the two cycles that 
each contain one of Pn and Pm or splits the single cycle 
that contains them both. This can be expressed formally 
as 



(Q\(n m)P) fl 
(Q\P) V2 



1 - 5 -" = fl if Pn Q ^Pm 
y j otherwise. 



(48) 



Here, the delta function takes care of the possibility that 
n = to, which corresponds to no transposition at all. 

Applying this result to Eqs. (EHJ (EH1> [with the right- 
hand side of Eq. I|48|) interpreted according to Q41JI] gives 







(e) 






fm 








(•) 




(o 




fo) 




(• 











FIG. 5: The operator XijX% ~ 1/4 — Si • Sj causes a re- 
arrangement of bonds, as shown in Fig. |31 Thus, its matrix 
element can be understood in terms of joining or splitting 
loops. 



(Qi44i p > 



(Q\P) 



3e tJ S a ' 



(49) 



where and ctj are unique labels for the loops passing 
through sites i and j. We have introduced the notation 



— 1 if i,j arc in different sublattices 
+ 1 otherwise. 



(50) 



To clarify, let us rederive Eq. I|49|l by considering ex- 
plicitly the effect of XijX% acting to the right on |P) 
and how the overlap of the resulting state with (Q\ dif- 
fers from the loop structure of (Q\P). We distinguish 
between the AB case, in which the sites i and j are in 
different sublattices, and the AA/BB case, in which they 
are in the same sublattice. For AB operators, the pos- 
sibilities are illustrated in Fig. OJa-c). If the two sites 
belong to the same loop then either an offdiagonal oper- 
ation splits the loop, giving a contribution (1/2) (2 1 ) = 1 
[following Eq. I|47|l ] or a diagonal operation leaves the 
loop unchanged, also giving (1)(2°) = 1. If the two sites 
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belong to different loops then the offdiagonal operation 
joins the two loops, giving (1/2) (2 -1 ) = 1/4. Combining 
these results gives 



{Q\x%x%\p\ 



(Q\P) 



5°»<*j + _M 

4 



-6 ai 
4 



f 

4' 



(51) 



For AA/BB operators, the possible bond reconfigurations 
are illustrated in Figs. OJd) andOJe). The two sites are 
either in the same loop, giving (f/2)(2° — 2 1 ) = —1/2, 
or in different loops, giving (1/2) (2° - 2" 1 ) = 1/4. The 
result, 



(Q\x%x%\P) 
(Q\P) 



= ~-6 a * 
4 



4 
1 

4' 



(52) 



differs from Eq. (|51|l by only a sign. Using ey to account 
for the difference in sign yields Eq. (|49|l . 

The same formal manipulations that lead from Eq. H48|) 
to Eq. (|4^|l can be chained together to evaluate a long op- 
erator sequence one transposition at a time. For example, 
taking i = i n £ A and j = j m e B, we get 



{Q\OxVjX%\P) = (Q\Q\(Pn m)P) m 1 "*"- 

(Q\P) (Q\p) W 



(Q\o\n (Q\n 
(Q\n (Q\p) 



l-8p 



(53) 



(Q\d\p') i, 

(Q\P>) 4 1 



where P' = (Pn m)P describes the new bond configu- 
ration after one transposition. For arbitrary site indices, 
the expression reads 



{Q\6x? 3 x%\P) 
(Q\P) 



i,,, AQ\o\P) 

4 U + £yJ (Q\P) 



(Q\o\n i 

(Q\P>) 4 



1 + 35° 



(54) 



VI. 



SPIN CORRELATION FUNCTIONS 



So long as the ground state of the spin system is a 
global singlet, its wavefunction can be written as a linear 
superposition of valence bond states: \tp) — J^p cp\P). 
Accordingly, operator expectation values have the form 



(0) = 



(jW) _ Ep,qW(P,Q) 

WW 



(Q\o\p) 
(Q\P) 



Z PiQ W(P,Q) 



(55) 



with P,Qg iS/v and O some operator of interest. The 
quantity W{P, Q) = cqcp(Q\P) may be determined vari- 
ational!}' 11 or it may arise as a sampling weight in the 



context of a Monte Carlo projection schemed The ma- 
trix element (Q\0\P) / (Q\P) is related to the properties 
of the loops that are formed when the singlet tilings of 
the two states are superimposed. In fact, since its value 
depends only on the properties of the loops, the individ- 
ual valence bond states can be abstracted away entirely, 
leaving what is essentially a loop estimator for O; we de- 
note this function Oc- In this way of thinking, Eq. I|55|) 
can best be understood as (O) = (Oc)w, an ensemble 
average of Oc in the gas of fluctuating loops. 38 ' 39 

In this section, we want to determine the loop estima- 
tors for the spin operators • Sj and (Si • Sj)(Sk -Si). All 
the work to compute the second order result has already 
been done. Since • Sj = j — XuXij> comparison with 
Eq. 149(1 immediately yields 



(56) 



For the fourth order result, we begin by specializing 
Eq. l(S^|) to the case O — X°kiXki> which gives 

(QIX ^|^ X "-' |jP) = ±(l-e ij 36< a !)(l-e kl 36 a " a ') 
+ ^(l + £ij )e W (^-^). (57) 

Here, a labels the loops in QP and a' labels the loops in 
QP' where P' is the modified configuration after xtiXki 
has acted on |P). Making use of 



Of Of 
Xij Xij XklXkl 



i-ISi-Sj.-iSfc-S, 



(Si-SjXSfc-Sj) (58) 



and Eq. (|56() . we arrive at 

[(Si-s^Sfc-sOk 

r 3 , m , 



3 / i °i , , 

16 v ' 16 



(59) 



There is no contribution when 8 a ^ a i = S aiOIj = 0. Hence, 
the estimator is nonzero only if all four site indices belong 
to the same loop or if there are two indices in each of two 
loops. The possible configurations are shown in Fig. 
Note that if the vertices and (k, I) reside in different 
loops or if they are in the same loop but remain unlinked 
(i.e., there is a path along the loop from i to j that 
does not encounter k or I) then J Q » a j = S aiaj . This 
equality holds because the operation on (k, I) does not 
disrupt the loop structure at sites i and j. In this case, 
the bracketed term on the right-hand side is jgS aiaj . 

Otherwise, = 1 — 8 aiaj and the bracketed term 

is ^(1 -2S a ^). 

The undesirable aspect of Eq. (|59|l is that it is history- 
dependent. There are several cases to consider because a' 
references the loop configuration as it exists after applica- 
tion of the (fc, I) vertex. Ideally, we want to reexpress the 
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FIG. 6: The expectation value of e i:J e fc; (Si ■ Sj)(S fc ■ Si) 
is nonzero only if and (k, I) connect one or two loops. 

For one loop, the value is either 9/16 or —3/16 depending 
on whether k and / are in the (a) same or (b) different loop 
segments connecting i and j. For two loops, the value is either 
9/16 or 3/16 depending on whether the interactions leave the 
loops (c) disjoint or (d) connected. 



This follows because there are l? a ways to choose two 
sites from a loop of length L a . In the same way, 



M 4 



E 



Si -Sj) (Sfe -S;) can be computed using 
Eq. H60fl . The only complication is the case where i, j, k, I 
are all in the same loop a. Of the such configurations, 



Et 2 ^° 



1] 



(62) 



have weight —3/16. The counting reflects the fact that 
there are L a ways to fix i and, as j ranges over the loop, 
2j(L a — j) — 1 ways to place k and I in opposite loop 
segments. (This includes mixed cases of the form i = k 1 
j 7^ I which are A = 1 in nature, but omits those of the 
form i = k, j = I, which arc A = 0.) The remainder have 
weight 9/16. Thus, the total contribution is 



E 



9 /2 



3 \ 



2 >4 + 4)£^< (63) 



16 \ 3 Q 3 a J 16 \3 a 3 a 

3 9 
16 + 16 

which simplifies to 



15 
16 



E^« 



(64) 



estimator in a way that eliminates the a' labels. This 
can be accomplished — at the expense of introducing a 
new quantity A — as follows: 



[(s i .s j )(s k .s l )] i 



{l + 2A)8 a ^5 a ' ak 5 akai 



16 

_|_ _l^i a k^ajai _|_ fiaiatga-ja^ 

16 



(60) 



A = 0, 1 is a topological term whose value is nonzero 
when traversal of the loop reveals an antisymmetric per- 
mulation of the site labels i,j,k,l. It distinguishes the 
Fig. EJb) configuration (A — 1) from that of Fig. Efa) 
(A = 0). For spin correlations of order four and above, 
the loop estimator is no longer simply a function of the 
loop labels a t . 

We can make use of Eqs. I|56f) and l|60l) to calculate 
powers of M = EieA ^ — YljeB which on a bipartite 
lattice corresponds to the staggered magnetization. At 
second order, we have M 2 
Eq. (H, 



Ylij £ySi • Sj and hence, via 



Ml 



(Q\M 2 \P) = 3 
(Q\P) " 4 



N 

E £ c 



(61) 



VII. CUMULANT GENERATING FUNCTION 

Calculating spin correlation functions as we did in 
the previous section involves keeping track of all pos- 
sible ways that a given number of site index pairs, 
(i, j), (k, I), . . . , can be assigned to a set of valence bond 
loops. For each such assignment, we then have to deter- 
mine the contribution to the correlation function based 
on how the loops are reconfigured by the valence bond op- 
erators Xij Xijj XfciXfci' ■ ■ ■ • This brute force approach be- 
comes increasingly cumbersome (and errorprone) as the 
order of the correlation function becomes large. At sec- 
ond and fourth order, the problem is manageable. For 
S,; • S j , there are are two distinct configurations — one that 
contributes (when i and j are in the same loop) and one 
that does not (when i and j are in different loops) — and 
for (Si ■ Sj)(Sj - Si), there are eight configurations — four 
that contribute (those shown in Fig. |SJ| and four that do 
not. But at sixth order, there are already 33 configura- 
tions to consider — an intractable nightmare. 

In this section, we present a vastly more simple ap- 
proach that emerges from a deeper understanding of the 
relationship between correlation functions and loops. To 
start, we introduce a new operator 



Hi 



1 

4 



Si 



1 



(1 



ot^o 



(65) 
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)(4-l)-3 



We now introduce a generating function 



F[a] = log(Q|exp( ^ a i6 % J \P\ 



(72) 





)(3-2)-3 



1 

4 



FIG. 7: Two examples of the correlation function defined in 
Eq. I71|l acting at third order. 



which is designed to have a positive definite matrix ele- 
ment, irrespective of the sublattice membership of i and 
j. What motivates this definition is a desire to compen- 
sate for the asymmetry between Eq. (I36[) and Eqs. I|37|) 
and (|38[) . which describe the effect of Xij X% 011 a valence 
bond state. In the AB case (e^ = —1), the operator 
7*j = Xij Xij produces only a rearrangement of bonds 



and leaves the phase of the state unchanged, which is 
the desired state of affairs; in the AA/BB case = 1), 
the operator 7^ — \ — X^jXij 1S engineered to behave in 
exactly the same way, by cancelling the first term and 
changing the sign of the second term on the right-hand 
side of Eqs. I|37|l and lj38|l . Pictorially, this amounts to 
removing the identity part and changing the sign of the 
transposition part from — i to i on the right-hand side 
of the third (bottom) update rule in Fig. |3| The rear- 
rangement of individual singlets is governed by 



lii[i,3\ = [i,3\ 



lij[i,l][k,j] =7ife[i,i][fc,i] = ^[i,j][k,l] 
for sites i, k £ A and j,l € B. Hence, 



7Wm |P)=2 5 — - l \(Pnm)P) 1 
\P) = 2 & "^- 1 \{Pn Pm)P), 
\P) = 2 5 "-'"- 1 |(n m)P). 



(66) 
(67) 



(68) 
(69) 
(70) 



Since %j involves only the transposition of bond pairs, 
its effect on loops is limited to splitting one loop into two 
(when i and j belong to the same loop) and joining two 
loops into one (when i and j belong to different loops). 
Its correlation function obeys the simple rule 




(Q\%j k i---\P) 
(Q\P) 



= 2 AN °- n > 0, (71) 



where ANq represents the change in loop number after 
the loops have been split or joined at {k,l), etc. 

Two examples are given in Fig. \7\ 



Its diagramatic expansion in powers of the coupling a%j , 
shown in Fig. [HI is completely analogous to the Gold- 
stone diagrams familiar from standard many-body the- 
ory. (There are only connected diagrams because of the 
linked-cluster theorem.) The n th order derivatives of F 
are the cumulants of the correlators defined in Eq. 17111 , 




d n F 



(73) 



a=0 



and are related to them by 



J %3 — KJ l3 , 
Cijki — Cijki — CijCki-, 

Cijklmn — Cijklmn ^V/'CfcZmn 
— CijCklCmn- 



(74) 
(75) 

(76) 



What is important about the cumulants is that their sub- 
tractive terms lead to perfect cancellation (e.g., Cijki = 
CijCki so that Cijki — 0) whenever the 7- vertices connect 
the loops into a network that cannot be disentangled by 
cutting it along two loop segments. The exact statement 
is that a cumulant vanishes unless its configuration is ir- 
reducible in the standard many-body diagram sense. See 
Fig. El 

The cumulants, in turn, are related to the physical spin 
correlation functions by 



4 

SijSki 



Sk I Sij ran &m n Sij k I 
ra n • 



Cij 


— S^ki + 


Cijki 


= Sijki — 


ijklmn 


— Sijklmn 



(77) 
(78) 

(79) 



Here, we have used the shorthand 



S, 



ijkl-- 



N £ fcr--) 



(Q\(S i -S j )(S k -S l )---\P) 

(QW) 



(80) 



Thus, computing spin correlation functions at any order 
is just a matter of determining the relevant irreducible 
diagrams, computing the corresponding cumulants, and 
then solving for the desired term in Eqs. I|77|l - H79(l . By 
this means, we recover the loop estimators Eq. (|56|l and 

Eq. dnoi - 

It is also possible to compute even powers of the stag- 
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FIG. 8: The diagrammatic expansion of the cumulant generating function F [from Eq. 1721 1 is shown to third order in the 
coupling a,ij. The black lines denote valence bond loops and the red lines j operators. 



gered magnetization by summing over all the site indices 
]TC<,=M| + A 2 , 



J2C m =M 4 £ -(M|) 2 , 

ijkl 

Cmimn =M e c - 3M 2 £ M A C 

ijklrnn 



2\3 



2(Mi) 



(81) 
(82) 
(83) 



At second order in M (first order in 7^), there are two 
irreducible diagrams: a one- loop configuration of weight 
1 and a two-loop configuration of weight 1/4. Since there 
are l? a ways to select two sites in loop a and L a Lp ways 
to select one site in loop a and one in loop /3, we have 



a ^ a 

52X 



N 1 



N 2 



(84) 



A glance at Eq. (|61|1 confirms that this method produces 
the correct result. 

At fourth order in M (second order in jij), there are 
again only two irreducible diagrams: a one- loop (cross) 
configuration of weight —3/4 and a two-loop configura- 
tion of weight 3/16. We have already worked out the 
counting in Eq. l|fj2^l . Hence, 



Mi = -- 



16 



Z 2 ^ 



(Ml) 



U 



(85) 



which is the result of Eq. I|64|) . 

At sixth order in M (third order in 7y), there are five 
irreducible diagrams. The diagrams of weight 3/32 and 
—3/32 both number Y^a^f3^^a^p anc ^ * nus cance l eacn 
other, so we only need to consider the other three: a 
one-loop (asterisk) configuration of weight 3/2, a one- 
loop (railway tie) configuration of weight 3/4, and a two- 
loop configuration of weight —3/16. Counting the num- 
ber of site index arrangements for these three cases re- 
quires some work, but it is no different in principle from 
what we did in Eq. I|62|) . We skip the details and simply 
report the result: 



16 
15 
4 
" 3 



Li 

Li 



32 
— j 
15 



(86) 



r 6 2 T, L2 « L2 M+Li-4) 



3M^M£ - 2(M^.) 3 . 



This simplifies to 



Z(^ 



2^ 



81 
64 



Z^ 



(■>') 



(87) 



In terms of the QP cycle lengths (L a 



2,k n 



the 
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+ 

3_ 
16 





F: 




+ 



16 



d- 
Z>- 

+ + ^ + (£) 



FIG. 9: At each order, the contributions to F can be classified 
according to their irreducibility. Nonirreducible diagrams are 
those that can be stitched together from self-energy parts of 
lower order in 7. The self-energy parts, shown as shaded arcs, 
are built by removing a single loop segment from the diagrams 
in Fig.|8] The values listed next to each diagram represent the 
C cumulant value associated with that configuration. C = 
for the composite terms. 



magnetization formulas are 



3 

a 

»# = £(-io£+4£)+i5(x;*»Y 

a ^ a ' 

M^ = £(64^-40fc^)+8l(^) 3 



(88) 
(89) 

(90) 



k a ), 



Recall that in Sect. El we used n k = J2a=i^(^ 
the cycle length distribution function, to derive an ex- 
pression for the overlap (Q\P). It is easy to see that the 
second moment of n k is related to the presence of long- 
range antiferromagnetic order in the system: a loop av- 
erage of Eq. dSHJ gives (M 2 ) = J2k 3k 2 (n k ) w . Similarly, 
the fourth and sixth powers of the staggered magneti- 



TABLE III: Asymptotic behaviour of the loop density and 
staggered magnetization in the thermodynamic limit as a 
function of the loop distribution tail. 



(n k )w ~ k p 


0((N )/N) 


0({M 2 )/N 2 ) 


< p < 1 


1 

N 


1 


p = l 


log N 
N 


1 


l<p< 2 


1 

N' 2 -P 


1 


p = 2 


1 

log N 


1 

log N 


2 <p< 3 


1 


1 

NP-' 2 


V > 3 


1 


1 

N 


short-ranged 


1 


1 

N 



zation have terms involving the cycle-length correlations 
{n k n k ,) w and {n k n kl n k „) w . 

The key feature is the tail of the cycle-length distribu- 
tion. If it is too weak, then the loop gas is characterized 
by a macroscopic number of small loops. In the oppo- 
site limit, the system has a vanishing number density of 
loops, although each loop contains a nonvanishing frac- 
tion of all the spins. These system-spanning loops are 
the foundation of the long-range order. As an exam- 
ple, suppose that (n k )w ~ k~ p , with the normalization 

set by the constraint N — X)fcLi k(n k )w- The standard 
large- N summation rules are 



N 



-L-N^'P + O(N-P) if p < 1 
logAT + 0(l) ifp = l 

C(p) + 0(N 1 -p) if p > 1, 



(91) 



where £(p) is the Reimann Zeta function. These can 
be used to derive the asymptotic behaviour listed in Ta- 
ble ^Q] Figure ITU1 depicts the N — > 00 behaviour of the 
loop number density and the average staggered moment: 



(No) 1 <\ , C(f>) M 9 ^ 



N 
(M 2 ) 



fc=i 



(92) 



Efc 2 KV-|^j^(2-p). (93) 



fc=i 



Clearly, long-range order exists for a range of exponents 
< p < 2. The phase transition at the critical value 
p c = 2 should be visible as a sharp step in the function 
(M 2 ) /(N ) 2 and in the Binder ratios^ (M 4 ) /(M 2 ) 2 and 
(M 6 )/(M 2 ) 3 . 

In the p — > 00 limit, the cycle- length distribution be- 
comes sharply peaked at the minimum length, (n k )w = 
N6k,i, which is the result for a fixed dimer configura- 
tion. More generally, a generic VBS state will have 
(n k )w ~ d(N() — k), where No is the size of the pla- 
quet (the basic unit of translational symmetry breaking) . 
When the distribution is short-ranged — having either an 
upper cutoff or an exponentially suppressed tail — there 
is no magnetic order. 
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FIG. 10: The squared staggered magnetization and the loop 
density are plotted as a function of the exponent p for a cycle- 
length distribution of the form {n,k)w ~ fc _p . The two quan- 
tities are complementary indicators of the long range order 
below p — 2. 



Here, /i = (/ii,//2, ■ ■ ■ ,MaO € is a vector of indices 
describing the singlet /triplet character of each bond. The 
set of such states is overcomplete with respect to the full 
Hilbert space (i N N\ > 2 2N ). The restriction to AB 
bonds is still meaningful in the extended basis since for 
every index pair fjf' , v" there is at least one set of indices 
fijV, \x' -,v' (in fact there are exactly four) such that the 
equation 



c 3Xlj 'X lk 







(95) 



has a nontrivial solution. This is the singlet/triplet gen- 
eralization of Eq. (|25|l . 

To understand how triplet states evolve under the 
Heisenberg hamiltonian, we must reproduce the analysis 
of Sec. II VI but now with Eqs. (|10fl and 1)221) specialized to 



[X°-,X?/]|vac)=^|vac> 



(96) 



and 



VIII. VALENCE BOND COVERINGS OF THE 
FULL HILBERT SPACE 

In previous sections, we discussed the properties of the 
S = valence bond states. We now turn our attention 
to the full Hilbert space and consider extended valence 
bond states of the form 



/V 



\P;fl) = Ux 



(94) 



[X% , X£jx« + ] Ivac) = \ T x ^ Xfe/|vac) . (97) 
As before, the operator —Hij has both a diagonal and 
offdiagonal part. If there is a pre-existing singlet bond 
across the active sites, the bond configuration is left un- 
changed. If there is a pre-existing triplet bond, however, 
the state is annihilated. Otherwise, the bonds are re- 
arranged and the singlet/triplet labels reassigned in a 
spin-conserving fashion. The update rule for —Hij with 
i E A and j G B, illustrated in Fig. i|llfl . is given by 



»njf> 



6^°\P;p) 



-1,0, Mn+l) 

I 



i PPn 



L) \ 



»Miv)) 



if Pm = n 
otherwise 



(98) 



There are some important differences from the purely 
singlet-bond case. There is now the possibility of anihi- 
lating a state by acting directly on a triplet bond, unlike 
in the singlet sector, where (HijHi'jiHii/j" ■ ■ -)\P) ^ 
always. Also, when bonds of different species interact, 
the weight T x ^ v is potentially negative: the r operators 
obey the identity 



^ Jt f 

T T T 



and hence 



(t° t 1 


r 2 r 3 \ 


T 1 T° 


T 3 -T 2 


r 2 -r 3 


T° T 1 


\r 3 r 2 


-r 1 r° / 







(99) 



(100) 



where e A " yA is the alternating symbol. 

Overlaps of extended valence bond states can be com- 
puted as before, only now we must use the most general 
form of the anticommutation relationship in Eq. (|22|l . We 
find that 



{Q\P\P;p) 



N 

2 N -N -Q Rq 



(101) 



Here, Nq is the number of cycles in QP and K a is a 
chirality factor whose value is given by the trace of the 
product of t""'t' 1 " in proper order around each cycle. 
The form of K a follows from the "boxcar" property 

^tr(---r vt r A )r A ^'"' = tr(- • •t^t'V'V'''). (102) 

A 
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IX. STAGGERED MAGNETIZATION 
COMPONENTS 
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FIG 11: Summary of the update rules for triplet bonds. 
The number of cross hatches on each bond indicates its triplet 
character. 



For example, the loop QP = (1 2 4 5 3) • • • , highlighted 
in Fig. has associated with it a factor 



1 



K a = -tr( 



.^lt T Ml r t / 2t r /J2 T ^4t T A'4 T ^5t r M5 r l'3t r A'3 



(103) 

Since K a takes the values or ±1, the overlap is no 
longer guaranteed to be positive definite, as it was in 
the case of only singlet bonds. If we do not mix triplet 
species, however, positive definiteness is preserved. In 
that case, K a = 1 if the loop contains an even number of 
triplets and K a = if the loop contains an odd number. 
More generally, if we move around a given loop in proper 
order, encountering triplets of species a, b, c, . . . along the 
way, then K a takes the value 



one triplet 

(-l) x +\i t \6 ab two triplets 

(_l)i+l/V*>c three triplets 

(-l)M(6 ab 6 cd + S ad S bc - S ac S bd ) four triplets 



(104) 



where \ counts the number of nonzero entries (triplets) 
in jl that belong to the loop in question. In other words, 
the overall sign depends on how many triplets arise from 
\P; jl) and how many from \Q;v). 



J 



The staggered magnetization operator has the special 
property that it can always be written so as to induce 
no rearrangement of bonds. This is simply a matter of 
grouping its terms to match the particular bond tiling of 
the state it is acting on. That is, for any P £ Sn, 



N 



M a 



£s?-£^ = £(^ 



CO \ 



(105) 



j£B 



Since the bracketed term behaves according to Eq. JT3J, 
we see that the effect of M a on an S = valence bond 
state is to create a superposition of states, each of which 
has one bond promoted to an a-triplet: 



N 



M a \P)= l J2\P;ar l 



(106) 



n=l 



Here, \P;a n ) denotes the state \P) jl) with jl having all 
zero entries except for a in the n th slot. Since lone triplets 
always produce zero-contribution loops, we find that 



N 



m= (Q\M a \P) =iTj 



(Q\p) 



, > = o, (107) 



\ (Q\P) 



which is the expected result for rotationally invariant 
states \P) and \Q). Two triplets of the same species, 
however, if they originate in different states, give unit 
weight whenever they lie in the same loop. Thus, 



JV N 



(M«M b ) c = J2i: {Q;a { ^ K) = 5*» £ e (108) 



m— 1 n— 1 



The final equality holds because, in a cycle of length k a , 
there are k a ways to place the m th Q link and k a ways 
to place the n th P link. 

Correlation functions at fourth order can be computed 
in much the same way. Two applications of M a to a 
valence bond state yields 



(M a ) 2 \P) =N\P)- £|P;<w 
and the overlap of two such states is 



(109) 



(Q\(M a ) 2 (M b ) 2 \P) 
(Q\P) 



f(Q; \P) , (Q\P- 
I W) (Q\p) 



E E 



(110) 
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The bracketed terms in Eq. (|110fl each take the value 
— 1 if the two triplets are in the same loop and vanish 
otherwise. In the case of one triplet species (a = b), 
the final term contributes if there is either a single loop 
containing all four triplets or a pair of loops containing 
two triplets each. These configurations all have weight 
+1. Applying the appropriate counting arguments yields 

(M°)£ = N 2 + 2N^ k a (k a - 1) + J2 k l - l ? 

a a 

+ J2k a (k a -i)kp(kp-i) + 2j2k 2 a k 2 . (Ill) 

Note that the m ^ n and ml ^ n' constraints have no 
effect when m and m! (m and n') are in one cycle and 
n and n' (n and m!) are in another. Equation 11111) 
simplifies considerably to give 



(M 



a\4 
C 



(112) 



isotropy demands that ((M a ) 2 ) = -|(M 2 ), which is con- 
firmed by comparing Eqs. (|88|l and l|108(l . Second, the 
identity 

M 4 ^^{M a f + ^{M a f{M b f (117) 

requires that 3 x Eq. (|112fl and 6 x Eq. (|116[) sum to 
Eq. I|89|) . which is indeed true. In summary, the loops es- 
timators for components of the staggered magnetization 
at second and fourth order are 

(M a M b ) c = 5 ab k l, ( 118 ) 

a 

[{M«f{M b f] c ^E[( 1 -0^-(l + 20^)] 

a 

+ (1 + 2**) (r, kl 

(119) 



If there are two triplet species (a ^ b), the situation is 
slightly more complicated. When all four triplets are in 
the same loop, there are 

2k a X> - l)(n - 2) = hi - 2kl + h 2 a (113) 

n— 3 

configurations, having weight +1, in which the a and 
b triplets appear in consecutive order around the loop 
(aabb) and 

k a 

2k a J2(n-l)(k-n + l) = --kl + -k 4 a (114) 



configurations, having weight —1, in which they appear 
in alternating order (abab). [Note that Eqs. (|113[) and 
( 11141 sum to k^ika — l) 2 , as required.] There is also a 
positive contribution when the a and b triplets are paired 
(aa)(bb) in two different loops, but none when they are 
paired (ab)(ab). Hence, 



[{M a ) 2 {M b ) 2 ] c — N 2 + 2NY / k a (k a - 1) 

a 



J2 *a(*a - l)fc/j(*/J - 1), (H5) 



which simplifies to 
[(Ar) 2 (M b ) 2 ] £ = ^(^ 2 -^ 4 ) + (^fc 2 ) 2 . (116) 



It is important to verify that these expressions are 
compatible with those of Sects. IVII and IVIII First, spin 



X. THE NEEL STATE 

For a bipartite spin system, in which the A and B site 
labels designate genuine sublattices, the Neel state can be 
written in terms of any AB valence bond configuration. 
This freedom exists because the staggered moments can 
be organized into spin- up/spin-down pairs such that 

N 

\ r ) = n b h n h n i vac > = n b \A P j™) ^ 



j£B 



n=l 



for any PeS N . Since b^b] l = ^(x-j + itfj), the Neel 
state can be viewed as a superposition of states with a 
fixed (but arbitrary) bond configuration and total triplet 
number ranging from to N: 



N 



\R) 



2 N/2 



\P)+iJ2\P;3 n )- ^|P;3 m ,3, 



n=l 



T,n=l 



(121) 

which counts the num- 



If we define t a (p) 
ber of a-triplets in the index vector ft, then the overlap 
of \R) with an extended valence bond state IP; fx) can be 
written as 



(R\P;ft)=S tl<0 5 t2fi e i ^ 2 ^2- N ^ 



(122) 



The overlap vanishes if the state \P; fx) contains any 
1- or 2-triplets. Otherwise it is independent of P and 
constant up to an overall phase factor with an angle 
given by 7r/2 x the number of 3-triplets. The derivation 
of Eq. (|122|) relics on nothing more than the (species) 
orthogonality of bonds that connect the same two sites 
[expressed in Eq. (fTT)|) ] . There is no need to account for 
bond-configurational mismatch because the permutation 
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FIG. 12: Matrix elements of 7 taken between a valence bond 
state and the reference Neel state. Closed loops containing 
one triplet do not contribute, whereas open strings with one 
triplet contribute the same value as their singlet-bond-only 
counterpart. 



P in Eq. I|121|) can be chosen equal to the permutation 
in|P;/2). 

Matrix elements taken between the Neel state and a 
valence bond singlet state are particularly easy to com- 
pute. For example, applying Eq. I|68|) . we find that 



\P) 



(R\P) 



1-Sf 



(R\{Pn m)P) 



(123) 



where the ratio of overlaps on the right-hand side is 1 
since both (R\P) and (R\(Pn m)P) have the value 2~ N / 2 . 
In general, the matrix element picks up a factor of | 
for each offdiagonal 7 operation (keeping in mind that 
subsequent operations act on the reconfigured bonds). 
The evaluation rule can be expressed as 



n operators 

(r\iXt^-\p) 

(R\P) 



2 N o- 



(124) 



where N$ now refers to the number of closed loops that 
are formed by the bonds and 7 vertices. 



In analogy with Eq. I|101|) . the rules can be extended 
to cover states with arbitrary numbers of triplets: 



n operators 



(R\ ---IP; g) _ N -n TT v 

(R\P) ~ LL Ka 



(125) 



In this case, the index a ranges over both the closed loops 
and the open strings (iV«_p in number) that are formed. 
A closed loop involving k < n operators and bonds num- 
bered 1,2, .... A; has associated with it a factor 



Ka = itr(T ( V il T ( V ia 



(126) 



An open string involving k < n operators and bonds 
numbered 1, 2, . . . , k + 1 has a factor 

K a = ^tr^V'^V 12 •••t°V i *+ 1 ). (127) 

Note that any lone triplet in a closed loop gives K a = 0, 
whereas a lone 3-triplct in an open string gives K a = i. 
See Fig. El 



XI. SINGLET-TRIPLET GAP 

Consider a quantum spin system whose hamiltonian H 
has a global singlet ground state \ipo}- Starting from an 
arbitrary singlet trial state 



|v ,w ) = £ c w |p)j 



(128) 



the true ground state wavefunction can be obtained by 
projection: 



lim e 



-tH 1 /.trial 



\w 



(129) 



From the eigen-equation H\ipo) = Eo\ipo), the ground 
state energy can be isolated by acting from the left with 
an appropriate reference state: 



En = 



(R\4>o) 



(130) 



For our purposes, it is most convenient to use a Neel 
reference state \R), since it has a constant-magnitude 
overlap with all states in the full Hilbert space. Thus, 



E 



EpW(P) 



(131) 



with W(P) = c P (R\P) = 2- N / 2 c P . This formulation is 
similar to that of Eq. (|55|) . except that here the basic 
objects of interest are bonds rather than loops. We can 
view the ground state energy E = (Hg)w as the expec- 
tation value of the bond estimator for the hamiltonian, 
Hb, in a fluctuating gas of valence bonds. 
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Similarly, by constructing a trial wavefunction in the 
triplet sector, 



T — OO 



T = 



ivi rial ) = M 3 i^ rial ) = i E c ^ ial i p ; 3 «)> ( 132 ) 



n=l P 

we can project onto the lowest-energy triplet state: 



l^i 



lim e~ TH M z \^ 



trial \ 



(133) 



The structure of this projection is closely related to that 
of Eq. Ijl29(l . In the singlet case, the evolution operator 
consists of a long string of —Hij operators that succes- 
sively reconfigure the bond configuration \P) of the trial 
state. Here, we simply need to compute the result of 
the same strings operating on |P;3 n ). We know, how- 
ever, that the update rules for singlet and triplet bonds 
are identical except that a direct diagonal operation kills 
a triplet. It is straightforward to reinterpret the opera- 
tor string by designating one bond of the trial state as a 
triplet and tracing its evolution. As shown in Fig. ^3 the 
triplet will either run the gauntlet or die in the attempt. 
Accordingly, we can write 



JV 



(134) 



n=l P 



where g n (P) counts the number of surviving triplets 
whose final destination is the n th bond of configuration 
P. By definition, 



JV 



0<J29n(P)<N. 



(135) 



The energy of the low-lying triplet state is 



{R\H\^) E P W(P)j:Li9n(P)^^ 



E P W(P)Zn=l9n(P) 



Here we have made use of the fact that (R\P;3 r , 
i(R\P). Alternatively, E\ can be expressed as 



(136) 



Ex 



(137) 



ir 



since (R\H\P; 3 n )/(R\P; 3 n ) and (R\H\P) / (R\P) differ 
only when g„ = 0. Combining Eqs. (|131fl and (|137|l gives 



El — Eq = 



{Y, n 9nH B ) w - (j2n9n) w (H B ) 



(En 9n) 



(138) 



A numerical measurement of the singlet-triplet gap via 
Eq. (|138|) has the advantage that it can be carried out 
during a simulation of the S — ground state properties. 
This gives a substantial error cancellation in stochastic 




FIG. 13: (Top) The evolution operator (in series expansion) 
involves long series of — Hij operator strings that map VB con- 
figurations in the trial state to VB configurations in the true 
ground state through sequential reordering of the bonds. In 
the figure, the propagation index— increases from right to left 
along the horizontal axis. A set of singlet valence bonds are 
arranged along the vertical spatial axis. The red bars denote 
nearest neighbour interactions. (Middle) The same operator 
string can be reinterpreted as acting on a trial state with one 
triplet, marked here in purple. (Bottom) A triplet state that 
is acted on directly by a diagonal operation is annihilated. 



methods over techniques where Eq and E\ are computed 
individually and subtracted. In some systems, however, 
the lowest-energy singlet and triplet states may differ so 
substantially that the reweighting of the operator string 
becomes inefficient. In that case, so few triplets survive 
the projection that (E„ 9n) ~ 0. 



XII. SPIN STIFFNESS 

Again, let us imagine that the S = ground state is 
obtained via projection: 



= lim e~ r ' & |V' trial ). 



(139) 
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Now suppose that we introduce a twist field cj>(r) repre- 
senting a local rotation of the spins about the 3 axis. We 
will be concerned with the limit in which the field gra- 
dient V0(r) is small with respect to the lattice spacing. 
For concreteness, we endow the field with a single long- 
wavelength mode 4>{r) = 4>q + r • q. As a consequence, 



the relative spin rotation angles satisfy 9i 



■Tj)-q, 



Equation (|20|l then suggests that the hamiltonian — 
assuming it has the form of Eq. <|30(l — transforms as 



b /~tab 



-H -» -H[q] =-H + q a T a + -q a q b G 
where the gradient and Hessian are given by 



fa = \H & - ^ ) • e ° (44 + 4\m ) - ■ 



(140) 



(141) 



and 



G 



ab 



e a (r i -r j )-e b 



x (-xW+44)+- ( 142 ) 



Here, + • • • represents additional interaction terms in 
Kijki and beyond. Equation (|141f> is a spin current op- 
erator that changes the number of triplet bonds by ±1. 
(See Fig.O) It does so in proportion to the interaction 
strength and the projection of the bond length onto the 
axis of propagation. Equation l|142|) , on the other hand, 
changes the triplet count by or ±2. 



J 



Under the influence of the twist field, the ground state 
energy is 



E[q\ = lim 



(R\H[q\e- T "^ 



00 (R\e- TH W\ip ) 



(143) 



By virtue of the variational principle, E[q\ must be sta- 
tionary with respect to variation. Its convexity at q = 
is a measure of the spin stiffness. Derivatives d a = d/dq a 
can be handled using the identity 



dT 'e-^- T '^d a He- T '" = -d a {e- rfl ). 



(144) 



Hence, 
dE[q\ 



(R\f a \^) , (R\(H-E)\r) 



q=0 «) m) 

with \il> a ) = a |V>)| q=o given by 



= 0, (145) 



dr' e- ( - T - T ' )H f a e- T ' H \i; tli 



trial \ 



(146) 



The state \ip a ) evolves from the singlet trial state but 
differs from the projected ground state in that it is con- 
structed with an e a -directed triplet injected at all times 

t' < T. 

Making use of Eq. ^ 145(1 . we can show that the spin 
stiffness is given by 



where 



d 2 E[q\ 



dq a dq b 



{R\G ab \ii) (R\f a \t[j b ) {R\T b \ii a ) (R\(H - E)\ip ab 



q=0 



(R\rl>) 



\i> ab )= / dT' J dT" e -(T-T")Hfa e -(T"-T')Hfb e -T'H^tri^ 
'dT 1 fdT" e -(T-T")Hfb e -(T"-T')Hfa e -r'H^tri^ 



+ / dr' e-^- T ^ ft G ab e- r ' tt \^ tliaX ). 
Jo 



r 



(147) 



(148) 



The state \ip ab ) evolves from the singlet trial state 



trial \ 



case for the singlet-triplet gap measurement). The spin 



with two e a - and e -directed triplets injected at times r' stiffness is related to the fluctuations in energy of the fi- 

and t" . nal bond configuration arising from operator sequences 

t-i n rm , ^ i j. j ai • i in which neither triplet is annihilated. 

Equation (|147|) can be evaluated using the reweightmg ^ 

trick introduced in Sect. lXIl The same operator sequence 
used to generate the singlet ground state can be reinter- 
preted to include all possible triplet pairs (at all possible 
starting locations, not simply at the r = end as was the 
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FIG. 14: Several possible bond reconfigurations are shown, corresponding to action by each term in T a 
are suppressed. 



The relevant prefactors 



XIII. SUMMARY 

We have presented in detail a formal framework for 
organizing calculations in the valence bond basis. This 
approach is based on manipulations of valence bond cre- 
ation and annihilation operators, xfj an d Xiji that act 
between any two sites of the spin lattice. These operators 
are similar to those introduced in Ref. |3^ for fixed dimer 
configurations but are endowed with an expanded oper- 
ator algebra that is compatible with the overcomplete 
basis of all possible dimer configurations. 

We have focused on what we call the AB valence bond 
basis, in which half the lattice sites are assigned the la- 
bel A and the other half B and only dimers connecting 
sites with different labels are allowed. (The use of AB 
bonds only is routine in the literature, especially when 
the hamiltonian of interest does not contain explicitly 
frustrating AA/BB interactions.) The states in this re- 
stricted basis are uniquely characterized by permutations 
of the B site labels. Overlaps and matrix elements of such 
states are related in a systematic way to the cycle struc- 
ture of the permutations. We have shown, for example, 
that the presence of antiferromagnetism in a bipartite 
system is related to the long tail behaviour of the cycle 
length distribution (rik)w- 

Correlation functions of the isotropic spin interaction 
Si ■ Sj are related to the cumulants of the operator 



The various contributions 



% = k( 1 + e ij) - e ijX%X%- 
can thus be computed in a straightforward way from the 
set of connected Goldstone diagrams — albeit with the di- 
agrams interpreted quite differently than they are in the 
usual quantum many-body physics context. (Moreover, 
there is no need to worry about which of i and j are A or 
B sites. The different cases are all handled automatically 
by the sign factor e^ .) As an example of this approach, 
we have derived explicit expressions for second-, fourth-, 
and sixth-order expectation values of the staggered mag- 
netization. 



We have also emphasized that, when triplet bonds are 
fully accounted for, the valence bond basis spans the en- 
tire Hilbert space, not merely the S = subspace. The 
valence bond operators carry an index /i that specifies 
the singlet (fi — 0) or triplet (/i = 1, 2, 3) nature of each 
bond. We have extended the rules for computing over- 
laps and matrix elements to include states with triplet 
bonds. As an example of how this can be useful, we 
have derived expressions for the singet-triplet gap and 
the spin stiffness that can be interpreted a reweightings 
of the operator string in a projected singlet ground state. 
The reweighting formulation is ideally suited for use in 
valence bond projector Monte Carlo. 

This material is based upon work supported by the 
National Science Foundation under Grant No. DMR- 
0513930. 



APPENDIX A: 

A permutation P £ Sn is a bijective map that takes n 
to Pn (and Pn to n) for every 1 < n < N: 



P : 



(l Pl\ 

2 P2 

\N pn) 



(Al) 



An alternative notation for Eq. (|A1(I involves cycles of the 
form (1 2 3), which expresses the functional relationship 
1 — > 2, 2 — > 3, and 3 — > 1. Every permutation can be 
decomposed into a product of disjoint cycles 



#{P) 



(A2) 



a=l 



where #(-P) is the number of cycles in P, £ a is a repre- 
sentative number (the lowest, say) in a particular cycle, 
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and k a is its length. Two numbers n and m are in the 

same cycle of P iff there is some integer power k > such 

, p 
that P n = to, which we denote n ~ to. It follows that 

ri is in the cycle labelled a iff n ~ £ Q . 

The cycle decomposition given in Eq. (|A2I) is unique 
and requires the fewest possible number of cycles. One 
can, however, further decompose each term into several 
smaller cycles that are no longer disjoint (i.e., they may 
have elements in common). This process can always be 
carried further until the permutation is expressed entirely 
as a product of transpositions (cycles of length 2). Ac- 
cordingly, to understand how the cycle structure of QP' 
differs from that of QP, we need only consider the se- 
quence of changes induced by the chain of transpositions 
that separates P' from P. 

The result of a single transposition acting on P is 



/ 1 

2 



[n m)P = 



Pl\ 

P2 



Hence, 



Q(n m)P 



I'll m 



Pm n 



\ N PN/ 



(I QPl\ 
2 QP2 

Pn Qm 

Pm Qn 

V N PQNj 



(A3) 



(A4) 



As to how this differs from QP, there are three possi- 
bilities to consider. First, when n — to, the 2-cycle is 
just an identity element; the equalities QP = Q(n m)P 
and #(QP) = #(Q(n m)P) follow trivially. Otherwise, 
the outcome depends on whether Pn and Pm are in the 
same cycle. If they are, the transposition splits one cycle 
into two, 



QP = (1 8 3 Pn 2 7 4 Pm 6 5) • • • (A5) 
Q(n m)P = (1 8 3 Pm 6 5)(2 7 4 Pn) ■ ■ ■ , (A6) 



and #(Q(n m)P) - #(QP) = +1. If Pn and Pm are in 
different cycles, the transposition merges two cycles into 
one, 



QP = (1 3 7 Pn 4 6) (2 8 Pm 5) ■■ • (A7) 

Q(n m)P = (1 3 7 Pm 5 2 8 Pm 4 6) • • • , (A8) 

and m)P) — #(QP) = —1. These three cases can 

be summarized by 



#(Q(n m)P) - #(QP) + 5 n , m = 



or, equivalcntly, 



-1 if Pn ~ Pm 



— 1 otherwise 



(A9) 



(Q\(n m)P) 
(Q\P) 



G) 



_ J 1 if Pn Pto 
i otherwise. 



(A10) 
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